////////////////////////////////////////////////////////////////////////////////
// This program is free software; you can redistribute it                     //
// and/or modify it under the terms of the GNU General Public License         //
// version 2 (or higher) as published by the Free Software Foundation.        //
//                                                                            //
// This program is distributed in the hope that it will be useful, but        //
// WITHOUT ANY WARRANTY; without even the implied warranty of                 //
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU          //
// General Public License for more details.                                   //
//                                                                            //
// Written and (C) by                                                         //
// Engin Tola                                                                 //
//                                                                            //
// web   : http://cvlab.epfl.ch/~tola/                                        //
// email : engin.tola@epfl.ch                                                 //
//                                                                            //
// If you use this code for research purposes, please refer to the            //
// webpage above                                                              //
////////////////////////////////////////////////////////////////////////////////

#include <iostream>
#include "daisy.hpp"

using namespace std;
using namespace kutility;

daisy::daisy()
{
    m_verbosity = 2;
    m_image = 0;
    m_w = 0;
    m_h = 0;

    m_rad = 0;
    m_rad_q_no = 0;
    m_th_q_no  = 0;
    m_hist_th_q_no = 0;
    m_grid_point_number = 0;
    m_descriptor_size = 0;

    m_smoothed_gradient_layers = NULL;
    m_dense_descriptors   = NULL;
    m_grid_points = NULL;
    m_oriented_grid_points = NULL;

    m_scale_invariant = false;
    m_rotation_invariant = false;

    m_scale_map = NULL;
    m_orientation_map = NULL;
    m_orientation_resolution = 36;
    m_scale_map = NULL;

    m_level_sigmas = NULL;

    m_descriptor_memory = false;
    m_workspace_memory = false;
    m_descriptor_normalization_threshold = 0.154; // sift magical number

    m_disable_interpolation = false;

    m_nrm_type = NRM_PARTIAL;

    m_cube_size = 0;
    m_layer_size = 0;
}

daisy::~daisy()
{
    deallocate( m_image );

    if ( !m_workspace_memory ) deallocate( m_smoothed_gradient_layers );
    deallocate( m_grid_points, m_grid_point_number );
    deallocate( m_oriented_grid_points, g_grid_orientation_resolution );
    deallocate( m_orientation_map );
    deallocate( m_scale_map );
    deallocate( m_level_sigmas );
    if ( !m_descriptor_memory ) deallocate( m_dense_descriptors );
}

void daisy::set_parameters( double rad, int rad_q_no, int th_q_no, int hist_th_q_no )
{
    m_rad = rad;                   // radius of the descriptor at the initial scale
    m_rad_q_no = rad_q_no;         // how many pieces shall I divide the radial range ?
    m_th_q_no = th_q_no;           // how many pieces shall I divide the angular range  ?
    m_hist_th_q_no = hist_th_q_no; // how many pieces shall I divide the grad_hist
    m_grid_point_number = m_rad_q_no * m_th_q_no + 1; // +1 is for center pixel
    m_descriptor_size = m_grid_point_number * m_hist_th_q_no;

    for ( int i=0; i<360; i++ )
    {
        m_orientation_shift_table[i] = i/360.0 * m_hist_th_q_no;
    }
    m_layer_size = m_h*m_w;
    m_cube_size = m_layer_size*m_hist_th_q_no;
}

float* daisy::get_dense_descriptors()
{
    return m_dense_descriptors;
}

double** daisy::get_grid_points()
{
    return m_grid_points;
}

void daisy::reset()
{
    deallocate( m_image );
    deallocate( m_grid_points, m_grid_point_number );
    deallocate( m_orientation_map );
    deallocate( m_scale_map );
    if ( !m_descriptor_memory ) deallocate( m_dense_descriptors );
}

void daisy::release_auxilary()
{
    deallocate( m_image );
    deallocate( m_orientation_map );
    deallocate( m_scale_map );

    if ( !m_workspace_memory ) deallocate(m_smoothed_gradient_layers);

    deallocate( m_grid_points, m_grid_point_number );
    deallocate( m_oriented_grid_points, g_grid_orientation_resolution );
    deallocate( m_level_sigmas );
}

void daisy::compute_grid_points()
{
    double r_step = m_rad / m_rad_q_no;
    double t_step = 2*pi()/ m_th_q_no;

    if ( m_grid_points )
        deallocate( m_grid_points, m_grid_point_number );

    m_grid_points = allocate<double>(m_grid_point_number, 2);
    for ( int y=0; y<m_grid_point_number; y++ )
    {
        m_grid_points[y][0] = 0;
        m_grid_points[y][1] = 0;
    }

    for ( double r=0; r<m_rad_q_no; r++ )
    {
        int region = (int)(r*m_th_q_no)+1;
        for ( double t=0; t<m_th_q_no; t++ )
        {
            double y, x;
            polar2cartesian( (r+1)*r_step, t*t_step, y, x );
            m_grid_points[(int)(region+t)][0] = y;
            m_grid_points[(int)(region+t)][1] = x;
        }
    }

    if ( m_verbosity > 2 )
    {
        cout<<"[daisy] grid points:";
        display( m_grid_points, m_grid_point_number, 2 );
    }
}

/// Computes the descriptor by sampling convoluted orientation maps.
void daisy::compute_descriptors()
{
    if ( m_verbosity >= 0 )
    {
        cout<<"#######################################################################################\n";
        cout<<"as of version 1.5, compute_descriptors() does not return normalized descriptors.\n";
        cout<<"call normalize_descriptors() after compute_descriptors() if you want normalization.\n";
        cout<<"additionally, you may change the normalization algorithm by calling set_normalization()\n";
        cout<<"#######################################################################################\n";
    }
    if ( m_scale_invariant    ) compute_scales();
    if ( m_rotation_invariant ) compute_orientations();
    if ( !m_descriptor_memory ) m_dense_descriptors = allocate <float>(m_h*m_w*m_descriptor_size);

    memset(m_dense_descriptors, 0, sizeof(float)*m_h*m_w*m_descriptor_size);

    int progress_counter=0;
    progress_bar bar( 0, m_h, 20);
    bar.set_text("[daisy] ");

    int y, x, index, orientation;
#ifdef USE_OPENMP
#pragma omp parallel for private(y,x,index,orientation)
#endif
    for ( y=0; y<m_h; y++ )
    {
        if ( m_verbosity > 0 ) cout<<bar(++progress_counter)<<flush;
        for ( x=0; x<m_w; x++ )
        {
            index=y*m_w+x;
            orientation=0;
            if ( m_orientation_map ) orientation = m_orientation_map[index];
            if ( !( orientation >= 0 && orientation < g_grid_orientation_resolution ) ) orientation = 0;
            get_unnormalized_descriptor( y, x, orientation, &(m_dense_descriptors[index*m_descriptor_size]) );
        }
    }
}

void daisy::smooth_layers( float* layers, int h, int w, int layer_number, float sigma )
{
    int fsz = filter_size(sigma);
    float* filter = new float[fsz];
    gaussian_1d(filter, fsz, sigma, 0);
    int i;
    float* layer=0;
#ifdef USE_OPENMP
#pragma omp parallel for private(i, layer)
#endif
    for ( i=0; i<layer_number; i++ )
    {
        layer = layers + i*h*w;
        convolve_sym( layer, h, w, filter, fsz, layer );
    }
    deallocate(filter);
}

void daisy::save_descriptor( string filename, int y, int x, bool single_row)
{
    float* feat = &(m_dense_descriptors[(y*m_w+x)*m_descriptor_size]);
    if ( single_row ) save( filename, feat, 1, m_descriptor_size );
    else             save( filename, feat, m_grid_point_number , m_hist_th_q_no );
}

void daisy::save_descriptors_ascii( string filename )
{
    kutility::save(filename, m_dense_descriptors, m_h*m_w, m_descriptor_size );
}

void daisy::save_descriptors_binary( string filename )
{
    kutility::save_binary(filename, m_dense_descriptors, m_h*m_w, m_descriptor_size, 1, kutility::TYPE_FLOAT );
}

void daisy::normalize_partial( float* desc )
{
    float norm;
    for ( int h=0; h<m_grid_point_number; h++ )
    {
        norm =  l2norm( &(desc[h*m_hist_th_q_no]), m_hist_th_q_no );
        if ( norm != 0.0 ) divide( desc+h*m_hist_th_q_no, m_hist_th_q_no, norm);
    }
}
void daisy::normalize_full( float* desc )
{
    float norm =  l2norm( desc, m_descriptor_size );
    if ( norm != 0.0 ) divide(desc, m_descriptor_size, norm);
}
void daisy::normalize_sift_way( float* desc )
{
    bool changed = true;
    int iter = 0;
    float norm;
    int h;
    while ( changed && iter < MAX_NORMALIZATION_ITER )
    {
        iter++;
        changed = false;

        norm = l2norm( desc, m_descriptor_size );
        if ( norm > 1e-5 )
            divide( desc, m_descriptor_size, norm);

        for ( h=0; h<m_descriptor_size; h++ )
        {
            if ( desc[ h ] > m_descriptor_normalization_threshold )
            {
                desc[ h ] = m_descriptor_normalization_threshold;
                changed = true;
            }
        }
    }
}
void daisy::normalize_descriptors( int nrm_type )
{
    time_t st, en;
    if ( m_verbosity > 0 )
    {
        cout<<"[daisy] normalizing descriptors...";
        time(&st);
    }
    int number_of_descriptors =  m_h * m_w;
    int d;

#ifdef USE_OPENMP
#pragma omp parallel for private(d)
#endif
    for ( d=0; d<number_of_descriptors; d++ )
        normalize_descriptor( m_dense_descriptors+d*m_descriptor_size, nrm_type );
    if ( m_verbosity > 0 )
    {
        time(&en);
        cout<<" in "<<difftime(en,st)<<" seconds\n";
    }
}

void daisy::initialize_single_descriptor_mode()
{
    initialize();
    compute_smoothed_gradient_layers();
    compute_oriented_grid_points();
}

void daisy::initialize()
{
    if ( m_verbosity > 0 ) cout<<"\n[daisy] initializing...";

    int glsz = compute_workspace_memory();
    if ( !m_workspace_memory ) m_smoothed_gradient_layers = new float[glsz];

    float* gradient_layers = m_smoothed_gradient_layers;

    if ( m_verbosity > 3 )
    {
        cout<<"\n[initialize] saving input.bin\n";
        save_binary("input.bin",m_image, m_h, m_w, 1, TYPE_FLOAT);
    }

    layered_gradient( m_image, m_h, m_w, m_hist_th_q_no, gradient_layers );
    if ( m_verbosity == 4 )
    {
        for ( int ii=0; ii<m_hist_th_q_no; ii++ )
        {
            string file = "gradient_layers"+num2str(ii)+".bin";
            cout<<"[initialize] saving "<<file<<endl;
            save_binary(file,gradient_layers+ii*m_h*m_w, m_h, m_w, 1, TYPE_FLOAT);
        }
    }

    // assuming a 0.5 image smoothness, we pull this to 1.6 as in sift
    smooth_layers( gradient_layers, m_h, m_w, m_hist_th_q_no, sqrt(g_sigma_init*g_sigma_init-0.25) );

    if ( m_verbosity == 4 )
    {
        for ( int ii=0; ii<m_hist_th_q_no; ii++ )
        {
            string file = "sgradient_layers"+num2str(ii)+".bin";
            cout<<"[initialize] saving "<<file<<endl;
            save_binary(file,gradient_layers+ii*m_h*m_w, m_h, m_w, 1, TYPE_FLOAT);
        }
    }

    compute_level_sigmas();
    compute_grid_points();
    if ( m_verbosity > 0 ) cout<<" ok!\n";
}

void daisy::compute_level_sigmas()
{
    if ( m_level_sigmas == NULL )
    {
        // user didn't set the sigma's; set them from the descriptor parameters
        g_cube_number = m_rad_q_no;
        m_level_sigmas = allocate<double>(g_cube_number);

        g_grid_radius[0] = 0;
        double r_step = m_rad / m_rad_q_no;
        for ( int r=0; r< m_rad_q_no; r++ )
        {
            g_grid_radius[r] = r;
            m_level_sigmas[r] = (r+1)*r_step/2;
            if ( m_verbosity > 2 ) cout<<"sigma "<<r<<": "<<m_level_sigmas[r]<<endl;
        }
    }
}

void daisy::set_level_gaussians( double* sigma_array, int sz )
{
    g_cube_number = sz;

    if ( m_level_sigmas ) deallocate( m_level_sigmas );
    m_level_sigmas = allocate<double>(g_cube_number);

    for ( int r=0; r<g_cube_number; r++ )
    {
        m_level_sigmas[r] = sigma_array[r];
        if ( m_verbosity > 1 ) cout<<"[daisy] sigma"<<r<<": "<<sigma_array[r]<<endl;
    }

    for ( int r=0; r<m_rad_q_no; r++ )
    {
        double seed_sigma = (r+1)*m_rad/m_rad_q_no/2.0;
        g_grid_radius[r] = quantize_radius( seed_sigma );
        if ( m_verbosity > 1 )
        {
            cout<<"[daisy] seed : "<<seed_sigma<<endl;
            cout<<"[daisy] g_grid_radius["<<r<<"] = "<<g_grid_radius[r]<<" sigma: "<<m_level_sigmas[ g_grid_radius[r] ] <<endl;
        }
    }
}

void daisy::compute_histograms()
{
    time_t st, en;
    if ( m_verbosity > 0 )
    {
        cout<<"[daisy] compute_histograms ";
        time(&st);
    }

    int r, y, x, ind;
    float* hist=0;

    for ( r=0; r<g_cube_number; r++ )
    {
        float* dst = m_smoothed_gradient_layers+r*m_cube_size;
        float* src = m_smoothed_gradient_layers+(r+1)*m_cube_size;

#ifdef USE_OPENMP
#pragma omp parallel for private(y,x,ind,hist)
#endif
        for ( y=0; y<m_h; y++ )
        {
            for ( x=0; x<m_w; x++ )
            {
                ind = y*m_w+x;
                hist = dst+ind*m_hist_th_q_no;
                compute_histogram( src, y, x, hist );
            }
        }
    }
    if (m_verbosity > 0 )
    {
        time(&en);
        cout<<" in "<<difftime(en,st)<<" seconds\n";
    }
}

void daisy::compute_smoothed_gradient_layers()
{
    time_t st, en;
    if ( m_verbosity > 0 )
    {
        cout<<"[daisy] computing smoothed orientation layers...";
        time(&st);
    }

    float* prev_cube = m_smoothed_gradient_layers;
    float* cube = NULL;

    double sigma;
    for ( int r=0; r<g_cube_number; r++ )
    {
        cube = m_smoothed_gradient_layers + (r+1)*m_cube_size;

        // incremental smoothing
        if ( r == 0 ) sigma = m_level_sigmas[0];
        else         sigma = sqrt( m_level_sigmas[r]*m_level_sigmas[r] - m_level_sigmas[r-1]*m_level_sigmas[r-1] );

        int fsz = filter_size(sigma);
        float* filter = new float[fsz];
        gaussian_1d(filter, fsz, sigma, 0);

#ifdef USE_OPENMP
#pragma omp parallel for
#endif
        for ( int th=0; th<m_hist_th_q_no; th++ )
        {
            convolve_sym( prev_cube+th*m_layer_size, m_h, m_w, filter, fsz, cube+th*m_layer_size );
        }
        deallocate(filter);
        prev_cube = cube;
    }
    if ( m_verbosity > 0 )
    {
        time(&en);
        std::cout<<" in "<<difftime(en,st)<<" seconds\n";
    }

    if ( m_verbosity == 4 )
    {
        for ( int r=0; r<g_cube_number; r++ )
        {
            float* cube=m_smoothed_gradient_layers+(r+1)*m_cube_size;
            for ( int ii=0; ii<m_hist_th_q_no; ii++ )
            {
                string file = "cube"+num2str(r)+"_layer"+num2str(ii)+".bin";
                cout<<"[initialize] saving "<<file<<endl;
                save_binary(file,cube+ii*m_h*m_w, m_h, m_w, 1, TYPE_FLOAT);
            }
        }
    }
    compute_histograms();
}

void daisy::compute_oriented_grid_points()
{
    m_oriented_grid_points = allocate<double>(g_grid_orientation_resolution, m_grid_point_number*2 );

    for ( int i=0; i<g_grid_orientation_resolution; i++ )
    {
        double angle = -i*2.0*pi()/g_grid_orientation_resolution;

        double kos = cos( angle );
        double zin = sin( angle );

        double* point_list = m_oriented_grid_points[ i ];

        for ( int k=0; k<m_grid_point_number; k++ )
        {
            double y = m_grid_points[k][0];
            double x = m_grid_points[k][1];

            point_list[2*k+1] =  x*kos + y*zin; // x
            point_list[2*k  ] = -x*zin + y*kos; // y
        }
    }

    if ( m_verbosity > 1 )
    {
        // if want to display the descriptor structure
        int str_w = (int)(2*m_rad+1);
        int str_size = square( str_w );

        int* structure = allocate<int>( str_size );

        for ( int ori=0; ori<g_grid_orientation_resolution; ori++ )
        {
            memset( structure, 0, sizeof(int)*str_size );
            for ( int reg=0; reg<m_grid_point_number; reg++ )
            {
                int y = (int)(m_oriented_grid_points[ori] [2*reg  ] + m_rad);
                int x = (int)(m_oriented_grid_points[ori] [2*reg+1] + m_rad);

                structure[y*str_w+x] = reg+1;
            }
//          display( structure, str_w, str_w, 1, 1, 2, 2, "  " );
//          wait_key();
        }
        deallocate(structure);
    }
}

/// sets a custom grid
void daisy::set_grid_points()
{
    // I should implement this
    cout<<"[daisy] set_grid_points::this is not implemented yet \n";
    exit(1);
}

void daisy::smooth_histogram(float *hist, int hsz)
{
    int i;
    float prev, temp;

    prev = hist[hsz - 1];
    for (i = 0; i < hsz; i++)
    {
        temp = hist[i];
        hist[i] = (prev + hist[i] + hist[(i + 1 == hsz) ? 0 : i + 1]) / 3.0;
        prev = temp;
    }
}

float daisy::interpolate_peak(float left, float center, float right)
{
    if ( center < 0.0 )
    {
        left = -left;
        center = -center;
        right = -right;
    }
    assert(center >= left  &&  center >= right);

    float den = (left - 2.0 * center + right);

    if ( den == 0 ) return 0;
    else           return 0.5*(left -right)/den;
}

int daisy::filter_size( double sigma )
{
    int fsz = (int)(5*sigma);

    // kernel size must be odd
    if ( fsz%2 == 0 ) fsz++;

    // kernel size cannot be smaller than 3
    if ( fsz < 3 ) fsz = 3;

    return fsz;
}

void daisy::compute_scales()
{
    cout<<"###############################################################################\n";
    cout<<"# scale detection is work-in-progress! do not use it if you're not Engin Tola #\n";
    cout<<"###############################################################################\n\n";

    int imsz = m_w * m_h;

    if ( m_verbosity > 0 )
    {
        cout<<"[daisy] detecting scales...\n";
        cout<<"[daisy] k: "<<g_sigma_step<<" scale_st: "<<g_scale_st<<" scale_en:  "<<g_scale_en<<endl;
    }

    float sigma = pow( g_sigma_step, g_scale_st)*g_sigma_0;

    float* sim = blur_gaussian_2d<float,float>( m_image, m_h, m_w, sigma, filter_size(sigma), false);

    float* next_sim = NULL;

    float* max_dog = allocate<float>(imsz);

    m_scale_map = allocate<float>(imsz);

    memset( max_dog, 0, imsz*sizeof(float) );
    memset( m_scale_map, 0, imsz*sizeof(float) );

    int i;
    float sigma_prev;
    float sigma_new;
    float sigma_inc;

    sigma_prev = g_sigma_0;
    for ( i=0; i<g_scale_en; i++ )
    {
        sigma_new  = pow( g_sigma_step, g_scale_st+i  ) * g_sigma_0;
        sigma_inc  = sqrt( sigma_new*sigma_new - sigma_prev*sigma_prev );
        sigma_prev = sigma_new;

        if ( m_verbosity > 0 )
            cout <<"[daisy] [i = "<<i<<"/"<<g_scale_en<<"] smoothing: sigma_inc = "
                 <<sigma_inc<<", sigma_new = "<<sigma_new<<", fsz : "<<filter_size(sigma_inc)<<endl;

        next_sim = blur_gaussian_2d<float,float>( sim, m_h, m_w, sigma_inc, filter_size( sigma_inc ) , false);

#ifdef USE_OPENMP
#pragma omp parallel for
#endif
        for ( int p=0; p<imsz; p++ )
        {
            float dog = fabs( next_sim[p] - sim[p] );
            if ( dog > max_dog[p] )
            {
                max_dog[p] = dog;
                m_scale_map[p] = i;
            }
        }
        deallocate( sim );

        sim = next_sim;
    }

    blur_gaussian_2d<float,float>( m_scale_map, m_h, m_w, 10.0, filter_size(10), true);

#ifdef USE_OPENMP
#pragma omp parallel for
#endif
    for ( int q=0; q<imsz; q++ )
    {
        m_scale_map[q] = round( m_scale_map[q] );
    }

//    save( m_scale_map, m_h, m_w, "scales.dat");

    deallocate( sim );
    deallocate( max_dog );
}
void daisy::compute_orientations()
{
    cout<<"#####################################################################################\n";
    cout<<"# orientation detection is work-in-progress! do not use it if you're not Engin Tola #\n";
    cout<<"#####################################################################################\n\n";

    time_t sto, eno;
    time(&sto);

    if ( m_verbosity > 0 )
        cout<<"[daisy] starting orientation computation\n";

    assert( m_image != NULL );

    int data_size = m_w*m_h;
    float* rotation_layers = layered_gradient( m_image, m_h, m_w, m_orientation_resolution );

    m_orientation_map = zeros<int>( data_size );

    int ori, max_ind;
    int ind;
    float max_val;

    int next, prev;
    float peak, angle;

    int x, y, kk;

    float* hist=NULL;

    float sigma_inc;
    float sigma_prev = 0;
    float sigma_new;

    for ( int scale=0; scale<g_scale_en; scale++ )
    {
        sigma_new  = pow( g_sigma_step, scale  ) * m_rad/3.0;
        sigma_inc  = sqrt( sigma_new*sigma_new - sigma_prev*sigma_prev );
        sigma_prev = sigma_new;

        smooth_layers( rotation_layers, m_h, m_w, m_orientation_resolution, sigma_inc);

        for ( y=0; y<m_h; y ++ )
        {
            hist = allocate<float>(m_orientation_resolution);

            for ( x=0; x<m_w; x++ )
            {
                ind = y*m_w+x;

                if ( m_scale_invariant && m_scale_map[ ind ] != scale ) continue;

                for ( ori=0; ori<m_orientation_resolution; ori++ )
                {
                    hist[ ori ] = rotation_layers[ori*data_size+ind];
                }

                for ( kk=0; kk<6; kk++ )
                    smooth_histogram( hist, m_orientation_resolution );

                max_val = -1;
                max_ind =  0;
                for ( ori=0; ori<m_orientation_resolution; ori++ )
                {
                    if ( hist[ori] > max_val )
                    {
                        max_val = hist[ori];
                        max_ind = ori;
                    }
                }

                prev = max_ind-1;
                if ( prev < 0 )
                    prev += m_orientation_resolution;

                next = max_ind+1;
                if ( next >= m_orientation_resolution )
                    next -= m_orientation_resolution;

                peak = interpolate_peak(hist[prev], hist[max_ind], hist[next]);
                angle = (max_ind + peak)*360.0/m_orientation_resolution;

                int iangle = int(angle);

                if ( iangle <    0 ) iangle += 360;
                if ( iangle >= 360 ) iangle -= 360;


                if ( !(iangle >= 0.0 && iangle < 360.0) )
                {
                    angle = 0;
                }

                m_orientation_map[ ind ] = iangle;
            }
            deallocate(hist);
        }
    }
//    save( m_orientation_map, m_h, m_w, "orientations");

    deallocate( rotation_layers );
    time(&eno);

    if ( m_verbosity > 0 )
        cout<<"[daisy] finished orientation computation in "<<difftime(eno,sto)<<" seconds\n";

    compute_oriented_grid_points();
}

void daisy::set_descriptor_memory( float* descriptor, long int d_size )
{
    assert( m_descriptor_memory == false );
    assert( m_h*m_w != 0 );
    assert( d_size >= compute_descriptor_memory() );

    m_dense_descriptors = descriptor;
    m_descriptor_memory = true;
}
void daisy::set_workspace_memory( float* workspace, long int w_size )
{
    assert( m_workspace_memory == false );
    assert( m_h*m_w != 0 );
    assert( w_size >= compute_workspace_memory() );

    m_smoothed_gradient_layers = workspace;
    m_workspace_memory = true;
}



